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ABSTRACT 

Observations indicate that mass accretion rates onto low-mass protostars are gener- 
ally lower than the rates of infall to their disks; this suggests that much of the protostel- 
lar mass must be accreted during rare, short outbursts of rapid accretion. We explore 
when protostellar disk accretion is likely to be highly variable. While constant a disks 
can in principle adjust their accretion rates to match infall rates, protostellar disks are 
unlikely to have constant a. In particular we show that neither models with angular 
momentum transport due solely to the magnetorotational instability (MRI) nor gravi- 
tational instability (GI) are likely to transport disk mass at protostellar infall rates over 
the large range of radii needed to move infalling envelope material down to the central 
protostar. We show that the MRI and GI are likely to combine to produce outbursts of 
rapid accretion starting at a few AU. Our analysis is consistent with the time-dependent 
models of Armitage, Livio, & Pringle (2001) and agrees with our observational study 
of the outbursting object FU Ori. 

Subject headings: accretion disks, stars: formation, stars: pre-main sequence 



1. Introduction 

The standard model of low-mass star formation posits the free-fall collapse of a protostellar 
molecular cloud core to a protostar plus disk during times of a few times 10 5 yr (e.g., Shu, Adams, 
& Lizano 1987), consistent with the statistics of protostellar objects in Taurus (Kenyon et al. 
1990, 1994). To build up a star over these timescales requires a time-averaged infall rate of order 
2 x 10~ 6 — lO _5 M0yr _1 , rates typically used in calculations of protostellar properties at the end 
of accretion (Stahler 1988; Hart mann, Cassen, &: Kenyon 1997). The numerical simulations of 



dynamic star cluster formation by lBate et al.l ([2003) found that stars and brown dwarfs formed in 



burst lasting ~ 2x 10 4 years, implying infall rates of ~ 10 4 to 10 5 M0yr . However, the accretion 
luminosity implied by such infall rates is considerably higher than typical observed protostellar 
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luminosities (Kenyon et al. 1990, 1994). This "luminosity problem" can be solved temporarily by 
piling up infalling matter in the circumstellar disk. Most of the mass must eventually be accreted 
onto the star, however. This requires major accretion events that are sufficiently short-lived that 
protostars are usually observed in quiescence. 

This picture of highly time-dependent accretion is supported by observations. Individual knots 
in jets and Herbig-Haro objects, thought to be the result of outflows driven by accretion energy, 
argue for substantial disk variability (e.g., Bally, Reipurth, & Davis 2007). The FU Ori objects 
provide direct evidence for short episodes of rapid accretion in early stages of stellar evolution, with 
accretion rates of 10 M Q yr _1 or more (Herbig 1977; Hartmann & Kenyon 1996), vastly larger 
than typical infall rates of < 1O _5 M0 yr _1 for low-mass objects (e.g., Kenyon, Calvet, & Hartmann 
1993; Furlan et al. 2008). 

The mechanism driving FU Ori outbursts is not yet clear. A variety of models have been 
proposed: thermal instability (TI; Lin & Papaloizou 1985; Bell &: Lin 1994); gravitational instability 
(GI; Vorobyov & Basu 2005, 2005); gravitational instability and activation of the magnetorotational 
instability (MRI; Armitage, Livio, & Pringle 2001; also Gammie 1999 and Book &: Hartmann 2005); 
and even models in which planets act as a dam limiting downstream accretion onto the star (Clarke 
h Syer 1996; Lodato & Clarke 2004). Our recent analysis based on Spitzer IRS data (Zhu et al. 
2007) led us to conclude that a pure TI model cannot work for FU Ori. 

In view of the complexity of the problem and the physical uncertainties we adopt a schematic 
approach. We start with the (optimistic) assumption that protostellar accretion can be steady. 
We then show that the GI is likely to dominate in the outer disk, while the MRI is likely to be 
important in the inner disk, and that mismatches between the GI and MRI result in non-steady 
accretion for expected protostellar infall rates. Q Our analysis agrees with the results found in 
the time-dependent outburst model of Armitage et al. (2001), and is consistent with our empirical 
analysis of the outbursting system FU Ori (Zhu et al. 2007). Although our results depend on 
simplified treatments of the GI and MRI, the overall picture is insensitive to parameter choices. 
We predict that above a critical infall rate protostellar disk accretion can be (relatively) steady; 
observational confirmation would help constrain mass transport rates by the GI and the MRI. 



2. Overview 

A disk with viscosity E| v will evolve at cylindrical radius r on a timescale 

tu ~ r 2 /u (1) 



1 GI and MRI here means turbulent states initiated by gravitational or magnetic instabilities, respectively. 

2 We use "viscosity" as shorthand for internal, localized transport of angular momentum by turbulence. We will 
also make the nontrivial assumption that external torques (e.g. MHD winds) can be neglected. 
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If this is comparable to the timescales over which mass is being added to the disk, then in principle 
the disk can adjust to an approximate steady state with infall balanced by accretion. To fix ideas, 
we assume that the disk beyond 1 AU is mostly heated by irradiation from the central protostar 
of mass M*, so that the temperature T oc r -1 / 2 . For a fully viscous disk, we adopt the usual 
parametrization of the viscosity v = ac 2 /f2, where c s is the sound speed (for a molecular gas) and 
Q is the (roughly Keplerian) angular velocity. Then 



rHl 

ac 2 e 



1.3 x 10 3 aZ\ M l T^l R AU yr , (2) 



where a_i = a/0.1 is the viscosity parameter M\ = M*/M & , Rau = r/AU, and T 30 o is the 
temperature at 1 AU in units of 300K. From this relation we see that a fully viscous disk might be 
able to keep up with mass infall over typical protostellar lifetimes of ~ 10 5 yr if the radius at which 
matter is being added satisfies Rau ^ 10 2 a_i. In a layered disk picture, the viscosity may need 
to be modified such that v = ac s h, where h should be the thickness of the active layer instead of 
the midplane scale height. However, this difference is significant only if the temperatures differ at 
the active layer and the midplane, which they do not unless there is some midplane viscosity. In 
any case, Eq. (2) prov ides an upper limit to Rau- Since typical observational estimates of infall 



radii are ~ 10- 100 AU (lKenyon et al.lll993l ). protostellar infall to a constant a disk is likely to pile 



up unless a is relatively large. 

A more serious problem is that protostellar disks are unlikely to have constant a. The best 
studied mechanism for angular momentum transport in disks, turbulence driven by the MRI (e.g., 
Balbus & Hawley 1998), requires a minimum ionization fraction to couple the magnetic fields to 
the mostly neutral disk. As substantial regions of protostellar disks will generally be too cold 
for thermal (collisional) ionization, ionization by nonthermal processes becomes important. This 
led Gammie (1996) to suggest a layered model in which non-thermally ionized surface layers are 
magnetically coupled while the disk midplane remain inert. We modify Gammie's analysis by 
assuming that the heating of the outer disk is not determined by local viscous dissipation but by 
irradiation from the central protostar, as above. 

The mass accretion rate in a layered disk is 

M = 67rrV 2 |-(2S a z,r 1/2 ) , (3) 

where T, a is the (one-sided) surface density of the active layers. Taking T, a = constant, and assuming 
that the disk temperature T oc r -1 / 2 , 

M = 5 x 10- 7 S 10 oT 3 ooa-i^C7 , (4) 

where £ioo = S a /100g cm -2 . 

Our nominal value of a = 0.1 may be reasonable for well-ionized regions, but it may be an 
overestimate for the outer regions of T Tauri disks (see §5.3). Also, the fiducial value for E a is based 
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upon Gammie's (1996) assumption of cosmic ray ionization, which may be an overestimate due to 
exclusion of cosmic rays by scattering and advection in the magnetized protostellar wind. X-rays 
provide a higher ionization rate near the surface of the disk but are attenuated more rapidly than 
cosmic rays (Glassgold & Igea 1999), yielding similar or smaller S a . Both calculations assume that 
absorption of ions and electrons by grains is unimportant, which is only true if small dust is highly 
depleted in the active layer (e.g., Sano et al. 2000, also Ilgner & Nelson 2006a,b,c). In summary, 
it is likely that the estimate in equation (j4]) is an upper limit, and thus it appears unlikely that 
the MRI can transport mass at r ~ a few AU at protostellar infall rates 2 x 10 -6 — 10~ 5 Mq yr _1 . 
MRI transport resulting from non-thermal ionization might however move material adequately in 
response to infall at r > 10 — 100 AU. 

On the other hand, if some nonmagnetic angular momentum transport mechanism can get 
matter in to r < 1 AU, thermal ionization can occur and activate the MRI. A minimum disk 
temperature is given by the effective temperature generated solely by local energy dissipation 

T > T e ff ~ 1600(MiM_ 5 ) 1/4 (R/0.2AU)~ 3/4 K , (5) 

where M_5 = M/1O _5 M0 yr -1 . Radiative trapping in an optically thick disk will make internal 
temperatures even higher. If T > 1400 K most of the silicate particles will evaporate, thus elimi- 
nating a major sink for current-carrying electrons. Therefore, high accretion rates can potentially 
activate the MRI on distance scales of order 1 AU or less. 

If magnetic angular momentum transport is weak then mass will accumulate in the disk until 
the disk becomes gravitationally unstable, at which point gravitational torques can transfer mass 
inward. GI alone may result cause accretion outbursts (Vorobyov & Basu 2006, 2008), although 
the details of disk cooling are crucial in determining if such bursts actually occur due to pure GI. 
Moreover the GI may be unable to drive accretion in the inner disk. GI sets in when the Toomre 
parameter 

where we have set the epicyclic frequency k ~ Q, appropriate for a near-Keplerian disk. At small 
radius O and c s will be large and therefore £ must also be large if we are to have GI. Since rapid 
accretion causes significant internal heating (compared to heating by protostellar irradiation), large 
surface densities imply significant radiative trapping, raising internal disk temperatures above the 
effective temperature estimate above. Thus, when considering rapid mass transfer by GI, either in 
a quasi-steady state or in bursts, it is necessary to consider thermal MRI activation in the inner 
disk. 

The above considerations suggest that the only way low-mass protostellar disks can accrete 
steadily during infall is if a smooth transition can be made from the GI operating on scales of 
~ 1 — 10 AU to the thermally-activated MRI at smaller radii. To test this idea, we have constructed 
a series of steady-state disk models with realistic opacities. We compute both MRI and GI steady 
models and then investigate whether a smooth, steady, or quasi-steady transition is likely. Our 
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results indicate that making the optimistic assumptions of steady GI and MRI accretion results in 
a contradiction for infall rates thought to be typical of low-mass protostars. 



3. Methods 

We compute steady disk models employing cylindrical coordinates (r, z) , treating radiative 
energy transport only in the vertical (z) direction. Energy conservation requires that 

where M* is the central star's mass and we have assumed that the disk is not so massive as to make 
its rotation significantly non-Kepler ian. Balance between heating by dissipation of turbulence and 
radiative cooling requires that 

where 

v = ac 2 jn (9) 

and 

dr = pndz , (10) 

and k is the Rosseland mean opacity. We have updated the fitting formulae provided by Bell & Lin 
(1994) for the Rosseland mean opacity to include more recent molecular opacities and an improved 
treatment of the pressure-dependence of dust sublimation (Zhu et al. 2007). The new fit and a 
comparison with the Bell & Lin (1994) opacity treatment is given in the Appendix. 



Convection has not been included in our treatment. iLin &: Papaloizoul ([1980 ) show that for a 
power law opacity (k = kqT 13 ), convection will occur when (5 > 1. Our opacity calculations show 
that 0>1 only occurs for T > 2000ET. As our steady-state analysis depends upon disk properties 
for T < 1400 K, the neglect of convection will not affect our results (see also Cassen 1993). 

We ignore irradiation of the disk by the central star, as we are assuming high accretion rates 
and a low central protostellar luminosity. The diffusion approximation (equation [8]) is adequate 
since the disk is optically thick at the high mass accretion rates (M > 10~ 7 M & / yr) we are interested 
in. 

We also require hydrostatic equilibrium perpendicular to the disk plane, 

dP z _ GM*pz 
dz r 3 

and use the ideal gas equation of state 



(11) 



P=-pT. (12) 
I* 
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Given a viscosity prescription, equations © - (|12p can be solved iteratively for the vertical structure 
of the disk at each radius, resulting in self-consistent values of the surface density S, and the 
temperature at the disk midplane T c . 

In detail, we use a shooting method based on a Runge-Kutta integrator rather than a relaxation 
method (e.g., D'Alessio et al. 1998) to solve the two-point boundary value problem. Given a and 
M at r, we fix z = Z{ and set T = T e jf and r = 2/3 (this is adequate in the absence of significant 
protostellar irradiation), then integrate toward the midplane. We stop when the total radiative 
flux = cT^fj at z = Zf. In general Zf ^ 0; we alter the initial conditions and iterate until Zf = 0. 

For an MRI active disk we fix a = olm, assuming the disk is active through the entire column. 
We then check to see if thermal ionization is sufficient or if the surface density is low enough that 
non-thermal ionization is plausible. The exact temperatures above which MRI activity can be 
sustained are somewhat uncertain; here we assume the transition occurs for a central temperature 
of 1400 K, when the dust grains that can absorb ions and electrons and thus inactivate the MRI 
(e.g., Sano et al. 2000) are evaporated. We set au = 10 -2 to 10 _1 to span a reasonable range 
given current estimates (see §5.3). 

For simplicity we neglect the possible presence of an actively accreting, non-thermally-ionized 
layer. This omission will not affect our results at high accretion rates, for which the layered con- 
tribution is unimportant (equation 2]) ; our approximation then breaks down for M < 1O _6 M0 yr _1 
for large values of S a and olm- 

For the steady GI disk models a is not fixed. Instead we start with a large value of a = olq 
and then vary olq until Q = 2. The adoption of the local treatment of GI energy dissipation 
requires some comment. Since gravity is a long-range force a local viscous description is not 
generally applicable (Balbus & Papaloizou 1999). However, as Gammie (2001) and Gammie & 
Johnson (2003) argue, a local treatment is adequate if A c = 2c 2 /(GS) = 2itHQ < r; here A c is 
the characteristic wavelength of the GI. More broadly, our main result involves order-of-magnitude 
arguments; that is, as long as inner disks must be quite massive to sustain GI transport, and as long 
as there is some local dissipation of energy as this transport and accretion occurs, steady accretion 
will not occur for a significant range of infall rates. To change our conclusions dramatically, one 
would need to show that the GI causes rapid accretion through the inner disk without substantial 
local heating. We return to this issue in §5.1. 



4. Results 

Figures la-d show steady disk results for a central star mass of 1M© and accretion rates of 10~ 4 , 
10 -5 , 10~ 6 , and 10~ 7 MQyr _1 . Proceeding counterclockwise from upper left, the panels show the 
central disk temperature, olq, the one-sided surface density E = dzp, and the viscous timescale 
r 2 /u as a function of radius. The solid curves show results for pure-GI models, while the dotted 
and dashed curves show results for a« = 0.1 and 0.01, respectively. 
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The upper left panels show that the central temperatures rise more dramatically toward small 
radius in the GI models than in the MRI models. The GI models have higher temperatures 
because their higher surface densities lead to stronger radiative trapping. The GI solutions in 
these high-temperature regimes are unrealistic because they assume the MRI is absent, when it 
seems likely the MRI will in fact be active. These high temperature states do, however, suggest 
the possibility of thermal instability in the inner disk at high accretion rates, especially as the 
solutions near ~ 3000 K represent unstable equilibria (e.g., Bell & Lin 1994; §5). We consider the 
MRI models to be inconsistent at T < 1400 K (collisional ionization would be absent) and when 
£ > T, a < lOOgcm" 2 . 

Can a smooth or steady transition between MRI and GI transport occur? The transition region 
would be the "plateau" in the temperature structure which occurs near T ~ 1400 K (see Figure 1). 
This plateau is a consequence of the thermostatic effects of dust opacity, which vanishes rapidly 
at slightly higher temperatures. A small increase in temperature past this critical temperature 
causes a large decrease in the disk opacity and thus the optical depth; this in turn reduces the 
radiative trapping and decreases the central temperature. Thus disk models tend to hover around 
the dust destruction temperature over roughly an order of magnitude in radius, with the plateau 
occurring farther out in the disk for larger accretion rates. Since the plateau is connected with the 
evaporation of dust, it corresponds to a region where we might expect MRI activity!! 

First consider the case M = 1O~ 4 M0 yr -1 (upper left corner of Figure 1). The plateau region is 
very similar in extent for all models. More importantly, oq ~ 10 -2 in this region, and so the surface 
densities of the GI and qm = 10~ 2 models are nearly the same. This suggests that a steady disk 
solution is plausible with a transition from GI to MRI at a few AU for these parameters. Depending 
on the precise thermal activation temperature for the MRI, a smooth transition at around 10 AU 
might also occur for au = 0.1. 

Next consider the case M = lO _5 M0yr _1 (upper right corner of Figure 1). Here oq ~ 10 -3 in 
the plateau region, with resulting surface densities much higher than for either of the MRI cases. 
This discrepancy in a and £ between the two solutions makes a steady disk unlikely. A small 
increase in surface density in a GI model near the transition region, resulting in increased heating 
and thus thermal activation of the MRI, would suddenly raise the effective transport rates by one or 
two orders of magnitude, depending upon au- The result would be an accretion outburst. This is 
qualitatively the same situation as proposed for outbursts in dwarf novae, where thermal instability 
is coupled to an increase in a from the initial low state to the high state (similar to what Bell & 
Lin 1994 adopted to obtain FU Ori outbursts). Our inference of non-steady accretion also agrees 
with the time-dependent one-dimensional models of Armitage et al. (2001) and of Gammie (1999) 



3 There will be hysteresis because the dust size spectrum in a parcel of gas will depend on the parcel's thermal 
history. Heating the parcel destroys the dust and the accumulated effects of grain growth. Cooling it again would 
presumably condense dust with small mean size (and therefore a strong damping effects on MHD turbulence). The 
opacity would then vary strongly with time as the grains grow again. These effects are not considered here. 
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and Book & Hartmann (2005), as discussed further in in §5. 

A similar situation holds at lO~ 6 M yr _1 , although the evolutionary (viscous) timescales of 
the GI model are of order 10 5 yr, comparable to protostellar infall timescales. At this infall rate, 
the disk would only amass ~ 0.1M Q = 0.1M*, and so the disk might not need to transfer this 
mass into the star to avoid GI. At 10~ 7 M Q yr~ 1 , evolutionary timescales become much longer 
than protostellar lifetimes, and become comparable to T Tauri lifetimes; disk material can pile up 
without generating GI transport and consequent thermal activation of the MRI. In addition, an 
ctjf = 0.1 value could result in a steady disk with surface densities low enough to be activated 
entirely by cosmic ray or X-ray ionization. This does not mean, however, that T Tauri disks do 
not have layered accretion, as the surface density distribution depends upon the history of mass 
transport. 

The results of our calculations are summarized in the M — r plane in Figure 2. The solid curves 
farthest to the lower right, labeled Rq, are the radii at which the pure Gl-driven disk would have 
a central temperature of 1400 K (at which temperature the dust starts to sublimate), and thus 
activate the MRI. Moving up and left, the solid curve labeled Rm denotes the radii at which a pure 
MRI disk of the given au would have a central temperature of 1400 K. When these two curves are 
close together, or cross, oq and au are similar, making possible a smooth transition between GI 
and MRI and thus steady accretion. In the (shaded) regions between these two curves the viscosity 
parameters diverge, making non-steady accretion likely. 

The radial regions at which we predict material will pile up, trigger the MRI, and result in 
rapid accretion lie in the shaded regions. The dotted curve shows Rq and Rm where the disk has a 
central temperature of 1800 K (at which temperature all dust has sublimated). Rq and Rm at 1800 
K are smaller than they are at 1400 K because of the plateau region discussed above. Thus if the 
MRI trigger temperature is higher the outbursts are expected to be shorter because the outburst 
drains the smaller inner disk (r < Rq) on the viscous timescale. 

Figure 2 indicates that non-steady accretion, with potential outbursts, is predicted to occur 
for infall rates < lO~ 5 M0yr _1 for a M = 0.01 and < 10~ 4 M o yr~ 1 for a M = 0.1. As described 
above, for M < 1O _6 M0 yr -1 outbursts are unlikely, simply because the transport timescales are 
too long. Outbursts are expected to be triggered at r ~ 1 — 10 AU for protostellar infall rates 
~ 10~ 5 — lO~ 6 M0yr _1 . These predictions are relatively insensitive to the precise temperature of 
MRI activation; the dotted curves in Figure 2 show the results for a critical MRI temperature of 
1400 K, which simply shift the regions of instability to slightly smaller radii without changing the 
qualitative results. 

The other shaded band in Figure 2 denotes the region where thermal instability might occur. 
The two limits correspond to the two limiting values of the "S curve" (e.g., Faulkner, Lin, & 
Papaloizou 1983) at which transitions up to the high (rapid accretion) state and the low (slow 
accretion) state occur. 

Figures 3-6 show results for central star masses of 0.3 and O.O5M , respectively. The predictions 
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are qualitatively similar to the case of the 1M Q protostar, with the exception that thermal instability 
is less likely for the brown dwarf. This also implies generally unstable protostellar accretion for 
more massive protostars during the time that they are increasing substantially in mass. 

Much of the overall behavior of our results derive from the general property that disk tempera- 
tures rise strongly toward smaller radii. For optically-thick viscous disks, the central temperatures 
are proportional to 

T c ~ T eff r^ oc MV4 r -3/4 (KRS) i/4 , (is) 

where r is the vertical optical depth. Thus, even changes in surface density for differing values of 
a result in modest changes in radii where a specific temperature is achieved. Changing the mass 
accretion rate has a bigger effect, because SocM. 

5. Discussion 

Our prediction of unsteady accretion during protostellar disk evolution is the result of the 
inefficiency of angular momentum transport of the two mechanisms considered here: the MRI, 
because of low ionization in the disk; and the GI, because it tends to be inefficient at small radii, 
where and c s will be large, forcing £ to be large. To provide a feeling of just how large the 
surface density must be for Q = 2 in the inner disk, at accretion rates of 10~ 4 and 10 _5 Mq yr" 1 for 
the 1M star the disk mass interior to 1 AU would have to be ~ 0.6M Q and ~ 0.5M Q , respectively 
(Fig. 7), which are implausible large. At some point the disk must accrete most of its mass into 
the star, forcing the inner disk temperatures to be very large and thermally activating the MRI, 
resulting in outburst of accretion. Here we consider whether the assumptions leading to this picture 
are reasonable, then discuss applications to outbursting systems. 

5.1. Outbursts? 

Our inference of cycles of outbursts of accretion - piling up of mass by GI transport, followed 
by thermal triggering of the MRI - was found in the models of Armitage et al. (2001), as well as in 
the calculations of Gammie (1999) and Book & Hartmann (2005). We have also found outbursting 
behavior in time-dependent two-dimensional disk models, to be reported in a subsequent paper 
(Zhu, Hartmann, & Gammie 2009). Here we compare our results with those of Armitage et al. . 

Figure 8 shows the results of our stability calculations for parameters and opacities adopted by 
Armitage et al. : a central star mass of 1M , oim = 0.01, and an assumed triggering temperature for 
the MRI of 800 K. Armitage et al. found steady accretion at an infall rate of M = 3 x 1O _6 M yr _1 
but outbursting behavior at 1.5 x 10" 6 M Q yr _1 . This is reasonably consistent with our calculations; 
Rm and Rq are close together at M = 3 x lO _6 M0yr _1 and cross near M = 1 x lO _5 M yr _1 , 
suggesting stable accretion somewhere in this range. Armitage et al. find that the MRI is triggered 
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at about 2 AU, whereas our analysis (for M ~ lO _6 M0yr~ 1 ) would suggest a triggering radius 
of about 3 AU. Our ability to reproduce the results of Armitage et al. is adequate, considering 
that steady models do not precisely reproduce the behavior of time-dependent models, and that 
the form of cxq used by Armitage et al. is somewhat different from ours, though it still retains the 
feature of non-negligible GI only for small Q. 

Our finding of non-steady accretion is the result of assuming no other s i gnific ant level of 
angular momentum transport that is not due to GI or thermal MRI. iTerqueml (|2008l ) has shown 
that steady accretion is possible for a layered disk accreting at M = 10 _8 MQ/yr if there is a non- 
zero (non-gravitational) viscosity in disk regions below the surface active layers. Simulations have 
indicated that active layers can have an effect on non-magn etically active regions below, producing 



a Reynolds stress promoting accretion in the lower regions ([Fleming Stone 



2008 



2003 



Turner k, Sano 



Ilgner Nelson! 120081 ) . We argue that this effect is unlikely to be important for the much 



higher accretion rates considered here, simply because the amount of mass transfer that needs to 
occur is much higher than what is sustainable by a non-thermally ionized surface layer. It seems 
implausible that a small amount of surface energy and turbulence generation can activate a very 
large amount of turbulence and energy dissipation in a much more massive region. 



5.2. Local vs. non-local GI transport 

We have adopted a local formalism for GI whereas it has non-local properties. Furthermore, 
we have adopted azimuthal symmetry in calculating the dissipation of energy whereas energy will 
be deposited in nonaxisymmetric spiral shocks. Neither of these assumptions is strictly correct. 

Boley et al. (2006) performed a careful analysis of the torques in a three-dimensional model of a 
self-gravitating disk, including radiative transfer. They found that the mass transfer was dominated 
by global modes, but could be consistent with a locally-defined a(r,t). This result did not hold 
near the inner and outer edge of their disk, although this is not surprising as these regions were 
characterized by Q > 2 and thus one would not expect the GI to be operating. Boley et al. were 
unable to address whether energy dissipation was localized. Nevertheless it is difficult to imagine 
that gravitational instability could avoid some heating in regions with Q ~ 1, and only relatively 
small amounts of heating are required to activate the MRI at small radii. 

The details of the disk temperature structure near 1 AU must be found by three dimensional 
simulations of the GI with realistic cooling. The analysis presented here suggests that pure GI in 
the absence of MRI tends to lead to very long transport times in the inner disk, as required by 
our low values of oq. This presents two potential technical problems for a numerical investigation: 
first, numerical viscosity must be smaller than olq to follow the evolution; and second, the disk 
must be followed over long, evolutionary timescales. It will be challenging to follow the GI near 1 
AU numerically. 
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5.3. What is a M ? 



The magnetic transport rate olm is constrained by both observations and theory. A recent 
review of the observational evidence by King, Pringle, & Livio (2007) argues that olm must be 
large, of the order 0.1-0.4, based in part on observations of dwarf novae and X-ray binaries where 
there is no question of gravitational instability. Our own analysis required a ~ 0.1 in FU Ori (Zhu 
et al. 2007). 



On the theoretical side the situation is murky. Early calculations (jHawlev et al.lll996l ) sug- 
gested that fo r "shearing box" mode l s wit h zero mean azimuthal and vertical field au — 0.01. 
Recent work (jFromang Papaloizoul 120071 ). however, shows that olm does not converge in the 
sense that au — * as the numerical resolution increases. 

But are the zero mean field models relevan t to astrophysical disk s? Global disk simulations 



( Hirose et al. 2004 : McKinney Sz Narayan 2007 : Beckwith et al. 20081 ). local disk simulations in 



which the mean field is allowed to evolve b ecause of the boundary conditions ([Brandenburg et al 



19951 ). and observations of the galactic disk (j Vailed 120041 ) all exhibit a "mean" azimuthal field when 
an average is taken over areas of > H 2 in the plane of the disk. This suggests that the zero mean 
field local models are a singular case, and that mean azimuthal field models are most relevant to 
real disks (strong vertical fields would appear to be easily removed from disks acc ording to the 
plausible phenomenological argument originally advanced by Ivan Ballegooijenl (|l989l )). 



So what do numerical simulations tell us about disks with mean azimuthal field? Recent 
work shows that in this case the outcome depe nds on the magnetic Prandtl number Ptm = u/rj 
(jFromang et alj|2007l : iLesur fc Longarettil 120071 ) [y = viscosity and rj = resistivity) and that au is 
a monotonically increasing function of Ptm- This intriguing result, and the fact that YSO disks 
have Ptm *C 1 throughout (although more dimensionless parameters are required to characterize 
YSO disks, where the Hall effect and ambipolar diffusion can also be important), might suggest 
that olm should be small. But the numerical evidence also shows that olm depends on v in the 
sense that the dependence on Ptm weakens as v decreases. In sum, the outcome is not known as v 
drops toward astr ophysically plausible values. Mean azimuthal field models with effective Ptm ~ 1 
Guan et al.1 (120081 ) are also not fully converged; they show that olm increases, albe it slightly, as the 
resolution is increased. For a mean field with plasma j3 = 400 iGuan et al.l ((2008) find oim = 0.03 
at their highest resolut ion. In disks with an initi al strong azimuthal magnetic field in equipartition 
with thermal pressure. I Johansen Levinl (|2008l ) find a = 0.1 resulting from a combination of the 
Parker instability and an MRI-driven dynamo. 

Very small olm would pose a problem for T Tauri accretion. In the layered disk model, Gammie 
estimated the accretion rate to be 



M 



2 x 10 



-8 



«M \ 2 

0.01/ 



lOOg cm 



M Q yr 



-i 



(14) 



where S a is the surface density of the layer which is non-thermally ionized. Thus with om ^ 10 3 
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it would be difficult to explain typical T Tauri accretion rates. 

On the other hand «m ~ 0.1 could cause the outer disks of T Tauri stars to expand to radii of 
1000 AU or more in 1 Myr (Hartmann et al. 1998). There is no particular reason why the om ~0.1 
that we estimated for the thermally-ionized inner disk region in FU Ori should be the same as the 
effective a in the outer disks of T Tauri stars, which cannot be thermally ionized. 

5.4. Protostellar accretion 

Our models predict that most low-mass protostars will be accreting more slowly than matter is 
falling onto their disks. This is consistent with observational results, as outlined in the Introduction. 
The results of Armitage et al. (2001) suggested that steady accretion might be possible at ~ 
3 x 1O~ 6 M yr _1 and above (for 1M ). We find a different result because we adopt a significantly 
higher temperature for thermal MRI activation, closer to that required for dust evaporation. This 
means that our MRI triggering occurs at smaller radii, where the GI is less effective. It does 
seem likely that higher activation temperatures than the 800 K adopted by Armitage et al. are 
more plausible. Even if thermal ionization in the absence of dust is sufficient at around 1000 K in 
statistical equilibrium, ionization rates are so low that equilibrium is unlikely (e.g. Desch 1998). 
We also note that Armitage et al. were unable to obtain the high accretion rates and short outburst 
durations characteristic of FU Ori objects, but Book & Hartmann (2005) were able to reproduce 
the FU Ori characteristics better with a higher MRI activation temperature. 

At infall rates > lO~ 4 M0yr _1 , our models predict (quasi-) steady accretion (also Armitage et 
al. 2001); but such high rates are not expected to last long, perhaps only during an initial rapid 
phase of infall (Foster & Chevalier 1993; Hartmann et al. 1994; Henriksen, Andre, & Bontemps 
1997). Testing this prediction may be difficult as relatively few objects will be caught in this phase 
and they will likely be heavily embedded. 

At lower infall rates, Gl-driven accretion timescales are longer than evolutionary times and/or 
layered MRI turbulence may produce sufficient mass transport. Thus, we would not expect out- 
bursts for Class II (T Tauri) stars. 

5.5. FU Ori outbursts 

In our radiative transfer modeling of the outbursting disk system FU Ori (Zhu et al. 2007), 
we found that to fit the Spitzer Space Telescope IRS spectrum the rapidly-accreting, hot inner 
disk must extend out to ~ 1 AU, inconsistent with a pure thermal instability model. In contrast, 
the results of this paper suggest thermal MRI triggering can occur at a few AU, in much better 
agreement with observation. 

Our recent analysis of the silicate emission features of FU Ori (Zhu et al. 2008) also suggests 
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that the disk becomes dominated by irradiation rather than internal heating at distances of > 1 AU, 
but this is consistent with the results of this paper, as irradiation from the central disk can dominate 
local viscous dissipation if the disk is sufficiently flared. 

We also found that the decay timescales of FU Ori suggest ocm ~ 1CP 1 ; large values of olm 
are more likely to lead to outbursting behavior. High inner disk accretion rates also make thermal 
instability more likely very close to the central star; the presence or absence of this instability may 
account for the difference in rise times seen in some FU Ori objects (Hartmann & Kenyon 1996). 



6. Conclusions 



Our study predicts that the disk accretion of low-mass protostars will generally be unsteady 
for typical infall rates. During the protostellar phase, GI is likely to dominate at radii beyond 1 
AU but not at smaller radii; in contrast, rapid accretion should drive thermal activation of the 
MRI in the inner disk. Because of the differing transport rates comparable to typical infall values 
results in high inner disk temperatures sufficient to trigger the MRI. This is a general conclusion, 
though if the external disk accretion is driven by GI, the radius at which the MRI can be triggered 
thermally is much larger, because of the high surface density needed to produce a low value of 
Q. Furthermore, Gl-driving in the inner disk results in a low value of olq, much lower than the 
expected aj.f, for a wide range of M. The feature of mass accumulation at low external a followed 
by a change to a high inner viscosity is similar to thermal instability models (and also Armitage 
et al. 2001). Thermal instabilities may also occur in the inner disk at very high accretion rates, 
enhancing the potential for non-steady protostellar accretion. 



A. Appendix: Rosseland mean opacity 



The Bell ii Lin (1994) Rosseland mean opacity fit has been widely used to study high temper- 
ature accretion disks (CV objects, FU Ori, et al. ) for more than a decade, with opacities generated 
almost two decades ago. Our understanding of opacity sources (especially dust and molecular line 



D'Alessio et al. 


1998, 


2001; 


Zhu et al. 


2007) 



We have generated Rosseland mean opacity assuming LTE for a w ide range of temperature and 
pressure during our study of FU Orionis objects (|Zhu et alJl2007l . 12008). The molecular, atomic, and 
i onized gas opacities have been calculated using the Opacity Distribution Fun ction (ODF) method 
((Castelli KuruczlbooJ : Isbordone et alJbood : ICastellill200a : Izhu et al]l2007h which is a statistical 
approach to handlin g line blanketing when millions of lines are present in a small wavelength range 



Kur ucz et alJll974T ) . The dust opacity was derived by the prescription in iD'Alessio et al.l (|200ll ) 



Zhu et al.ll2008l ). Our opacity has been used not only to study FU Orionis objec ts but also to fit the 



gas opacity for Herbig Ae star disks constrained by interferometric observations (jTannirkulam et al 



2008). The opacities are shown in Figure [9]). Compared with I Alexander &: Ferguson! (|1994T ) or Zhu 
et al. (2007,2008), the Bell & Lin opacity lacks water vapor and TiO opacity around 2000 K and 
has a lower dust sublimation temperature. 



We have made a piecewise power-law fit to the Zhu et al. (2007, 2008) opacity (analogous to 
the Bell & Lin fit) to enhance computational efficiency (Table HJ see also Figured]). This speedup 
has been useful in performing the calculations of this paper, and is essential for our forthcoming 
two-dimensional hydrodynamic simulations of FU Ori outbursts (Zhu, Hartmann, & Gammie 2009). 

We acknowledge useful conversations with Ken Rice and Dick Durisen. This work was sup- 
ported in part by NASA grant NNX08A139G, by the University of Michigan, and by a Sony 
Faculty Fellowship, a Richard and Margaret Romano Professorial Scholarship, and a University 
Scholar appointment to CG. 
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Fig. 1. — Steady-state disk calculations for four accretion rates - 10 4 , 10 5 , 10 6 , and 
1O~ 7 M yr 1 , assuming a central star of mass 1M . The solid curves show solutions for Gl-driven 
accretion, as described in the text. The dashed and dotted curves yield results for steady disk 
models with a constant a = 10 -2 and 10 _1 , respectively (see text) 
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Fig. 2. — Unstable regions in the r — M plane for a 1M central star. The shaded region in the 
lower right shows where the central temperature of steady GI models exceeds an assumed MRI 
trigger temperature of 1400 K. The dotted curves show Rm and Rq (the boundaries of the shaded 
region; see text for definition) for an MRI trigger temperature of 1800 K. The shaded region in the 
upper left shows the region subject to classical thermal instability. 




Fig. 3. — Same as in Figured] but for a central star mass of O.3M0. 
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Fig. 4. — Same as in Figured] but for a central star mass of O.O5M0. 
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Fig. 5. — Same as figure [2] for 0.3M Q central star. 




Fig. 6. — Same as Figure [2] for the 0.05M Q central star. 
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Fig. 7. — The mass of the disk integrated between radius R and the outer radius of 20 AU for 
steady-state Q=2 disks, at four accretion rates - 10~ 4 , 10 -5 , 10~ 6 , and lO _7 M0yr _1 . The central 
star mass is 1M Q . 
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Fig. 8. — Same as Figure [2] for the parameters of Armitage et al. (2001) (see text) 
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Fig. 9. — Rosseland mean opacities: the dotted lines show the Bell & Lin (1994) fit, the solid curves 
show the detailed opacity calculation of Zhu et al. (2007, 2008) (solid line), and the dashed lines 
show the simple fit to Zhu et al. opacities (Table 1). 
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Table 1: Fit to Zhu et al. (2007, 2008) opacity 



log 10 T log 10 k comments 

< 0.03 log 10 P + 3.12 0.738 log 10 T- 1.277 grain opacity 

< 0.0281 log 10 P + 3.19 -42.98 log 10 T + 1.312 log 10 P + 135.1 grain evaporation 

< 0.03 log 10 P + 3.28 4.063 log 10 T - 15.013 water vapor 

< 0.00832 log 10 P + 3.41 -18.48 log 10 T + 0.676 log 10 P + 58.93 

< 0.015 log 10 P + 3.7 2.905 log 10 T + 0.498 log 10 P - 13.995 molecular opacities 

< 0.04 log 10 P + 3.91 10.191og 10 T + 0.382 log 10 P- 40.936 H scattering 

< 0.28 log 10 P + 3.69 -3.36 log 10 T + 0.928 log 10 P + 12.026 bound-free,free-free 

else a —0.48 electron scattering 



"with two additional condition to set the boundary: if log 10 k <3.586 log 10 T -16.85 and log 10 T < 4, log 10 k = 3.586 
log 10 T-16.85; if log 10 T < 2.9, log 10 k=0.738 log 10 T -1.277 



